Autor/a

Santiago García Sánchez

Introducción

Se comparan distintos paquetes para graficar curvas de supervivencia con R.

Los conjuntos de datos que se emplean para los ejemplos son los siguientes:

  • veteran (incluido en el paquete survival): contrasta los resultados de dos regímenes de tratamiento para el cáncer de pulmón.
  • bladder (paquete survival): contiene datos sobre la recurrencia del cáncer de vejiga.
  • pharmacoSmoking (paquete asaur): datos provenientes de un ensayo aleatorizado que compara la eficacia de la triple terapia frente al uso del parche para dejar de fumar.
Código
library(survival)
library(asaur)

veteran <- survival::veteran
bladder <- survival::bladder
pharmacoSmoking <- asaur::pharmacoSmoking

Con las funciones Surv(), survfit() y coxph() del paquete survival se ajustan distintas curvas de supervivencia:

  • Para el dataset veteran, se estiman funciones de supervivencia según tratamiento utilizando el método de Kaplan-Meier.
  • Para el dataset bladder, se obtiene la función de supervivencia con el estimador de Fleming-Harrington. Esta estimación será utilizada posteriormente para construir la función hazard estimada por el método de Nelson-Aalen.
  • Con el dataset pharmacoSmoking, se ajusta un modelo de hazards proporcionales de Cox considerando las variables grp (grupo asignado), ageGroup4 (edad) y employment (tipo de empleo). Luego, se obtiene la curva de supervivencia estimada por el modelo.
Código
km_fit <- survfit(Surv(veteran$time, veteran$status == 1) ~ veteran$trt,
                  type = "kaplan-meier", conf.int = 0.95, conf.type = "log-log",
                  error = "greenwood")

na_fit <- survfit(Surv(bladder$stop, bladder$event == 1) ~ 1,
                  type = "fh", conf.type = "log-log")

cox <- coxph(Surv(ttr, relapse) ~ grp + ageGroup4 + employment,
             data = pharmacoSmoking)
cox_fit <- survfit(cox)

R base y survival

Se realizan gráficos con R base y survival:

Código
plot(km_fit,
     main = "Función de supervivencia estimada por el método de Kaplan-Meier",
     xlab = "Tiempo (días)", ylab = "Probabilidad de supervivencia estimada",
     cex.main = 0.9, lty = 1, lwd = 2, cex = 1, cex.lab = 0.95, mark = 1,
     conf.int = TRUE, col = c("steelblue", "pink2"))

Código
plot(na_fit$time, main = "Función Hazard estimada por Nelson-Aalen",
     y = na_fit$n.event / na_fit$n.risk,
     type = "p", pch = 20, xlab = "Tiempo (días)", ylab = "Tasa")
lines(lowess(na_fit$time, na_fit$n.event / na_fit$n.risk), col = "violet")

Código
plot(cox_fit, main = "Modelo de regresión de Cox",
     cex.main = 0.9, lty = 1, lwd = 2, cex = 1, cex.lab = 0.95, mark = 1,
     col = "red3",
     ylab = "Probabilidad de supervivencia estimada", xlab = "Tiempo (días)")

Paquetes ggplot2 y ggfortify

Se utiliza la función autoplot.survfit() del paquete ggfortify y varias funciones del paquete ggplot2:

Código
library(ggplot2)
library(ggfortify)

autoplot(km_fit, conf.int.alpha = 0.2, censor.shape = "x") +
  labs(x = "Tiempo (días)",
       y = "Probabilidad de supervivencia estimada",
       title = "Función de supervivencia estimada por el método de Kaplan-Meier",
       color = "Tratamiento", fill = "Tratamiento") +
  theme(title = element_text(size = 10.5),
        plot.title = element_text(size = 11.5, face = "bold"))

Código
autoplot(cox_fit, conf.int.alpha = 0.2, conf.int.linetype = "dashed",
         conf.int.colour = "red", surv.colour = "red3", surv.size = 1.2,
         censor = FALSE)  +
  labs(x = "Tiempo (días)",
       y = "Probabilidad de supervivencia estimada",
       title = "Modelo de regresión de Cox")  +
  theme(title = element_text(size = 10.5),
        plot.title = element_text(size = 11.5, face = "bold"))

Paquete GGally

Se utiliza la función ggsurv() del paquete GGally, que genera curvas de supervivencia automáticamente:

Código
library(GGally)

ggsurv(km_fit)  +
  labs(x = "Tiempo (días)", y = "Probabilidad de supervivencia estimada",
       title = "Función de supervivencia estimada por el método de Kaplan-Meier",
       color = "Tratamiento", linetype = "Tratamiento")

Código
ggsurv(cox_fit, surv.col = "red3")  +
  labs(x = "Tiempo (días)",
       y = "Probabilidad de supervivencia estimada",
       title = "Modelo de regresión de Cox")

Paquete survminer

El paquete survminer ofrece la función ggsurvplot() para dibujar curvas de supervivencia:

Código
library(survminer)

ggsurvplot(km_fit, data = veteran, risk.table = TRUE, conf.int = TRUE,
           break.time.by = 50, xlim = c(0, 600), surv.median.line = "hv",
           ggtheme = theme_minimal())

Código
ggsurvplot(cox_fit, data = pharmacoSmoking, palette = "red3", risk.table = TRUE,
           conf.int = TRUE, surv.median.line = "hv", break.time.by = 50,
           ggtheme = theme_minimal())

Paquete ggsurvfit

Se construyen curvas de supervivencia con tablas de riesgo recurriendo a funciones del paquete ggsurvfit (principalmente ggsurvfit()).

Código
library(ggsurvfit)

survfit2(Surv(time, status) ~ trt, data = veteran) %>%
  ggsurvfit() +
  add_confidence_interval() +
  add_risktable(risktable_group = "risktable_stats") +
  scale_ggsurvfit() +
  add_risktable_strata_symbol(symbol = "\U25CF", size = 10)

Código
coxph(Surv(ttr, relapse) ~ grp + ageGroup4 + employment, data = pharmacoSmoking) %>%
  survfit2() %>%
  ggsurvfit() +
  add_confidence_interval() +
  add_risktable() +
  scale_ggsurvfit()  +
  scale_y_continuous(limits = c(0, 1))

Paquete rms

El paquete rms incluye las funciones psm() y npsurv() para modelos de supervivencia paramétricos y no paramétricos respectivamente. Las curvas ajustadas pueden graficarse con survplot() (visualizaciones estáticas) o survplotp() (visualizaciones interactivas).

Código
library(rms)

km_fit_rms <- npsurv(Surv(time, status == 1) ~ trt, data = veteran)

survplot(km_fit_rms,
         col = c("steelblue", "red"),
         xlab = "Tiempo (días)",
         ylab = "Probabilidad de supervivencia estimada")

Código
survplotp(km_fit_rms,
          col = c("steelblue", "pink2"),
          xlab = "Tiempo (días)",
          ylab = "Probabilidad de supervivencia estimada")
Código
survplot(psm(Surv(bladder$stop, bladder$event == 1) ~ 1),
         what = "hazard", xlab = "Tiempo (días)")

Paquete jskm

La función jskm() del paquete del mismo nombre construye un gráfico Kaplan-Meier junto con una tabla de riesgo.

Código
library(jskm)

jskm(km_fit, table = TRUE, ci = TRUE,
     main = "Función de supervivencia estimada por el método de Kaplan-Meier",
     ystrataname = "Tratamiento", ystratalabs = c("1", "2"),
     xlabs = "Tiempo (días)", ylabs = "Probabilidad de supervivencia estimada",
     label.nrisk = "N.º en riesgo")

Paquete KMunicate

El paquete KMunicate crea gráficos Kaplan-Meier utilizando un estilo que sigue una serie de recomendaciones:

Código
library(KMunicate)

KMunicate(km_fit, seq(0, max(veteran$time), length.out = 10))